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Abstract 

In this paper, we review recent developments on the derivation and properties of 
macroscopic models of collective motion and self-organization. The starting point is 
a model of self-propelled particles interacting with its neighbors through alignment. 
We successively derive a mean- field model and its hydrodynamic limit. The resulting 
macroscopic model is the Self-Organized Hydrodynamics (SOH). We review the 
available existence results and known properties of the SOH model and discuss it 
in view of its possible extensions to other kinds of collective motion. 
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1 Collective dynamics and self-organization 

Many different kinds of interacting agent systems can be observed in nature, such as 
bird flocks, fish schools, insect swarms, etc. They provide fascinating examples of self- 
organizing systems which are able to produce large scale stable coherent structures: a 
typical example of such a structure is a social insect nest (such as a termite nest). The 
structure scale exceeds the insect typical size by several orders of magnitude, and no 
agent in the community has the cognitive capacity of planning it. Therefore, it emerges 
as a product of the local interaction between the agents (here their ability to manipulate 
mud bullets in response to chemical signals deposited by the other agents), without the 
intervention of a leader [15]. Self-organization is ubiquitous and can be observed in the 
inanimate world (see e.g. the formation of galaxies, crystals, tornadoes), as well as in the 
living world (see examples above) and the social worlds (e.g. in traffic, crowds, opinion 
formation, finance, etc.). It appears at so many different scales that one can wonder 
whether the fate of the universe is a journey towards ultimate disorder, as pictured by 
Boltzmann in the concept of entropy. The concept of evolution to disorder and self- 
organization are somehow contradictory. The biologist J. Monod realized that these two 
observations are not easily reconciled and tried to do so in his famous essay 'Chance and 
Necessity' [JH] . 

Another concept related to self-organizing systems is that of criticality. Indeed, such 
systems exhibit phase transitions between a disordered state and self-organized ones (see 
a review in [57]). The 'thermodynamic variable' which induces phase transitions is often 
related to the intensity of the noise undergone by the agents in their motion. This provides 
an easy analogy to the temperature in classical thermodynamics. However, most often, 
the density or size of the system is another variable inducing phase transitions and more 
surprisingly, order appears at large densities. This is somewhat in contradiction with 
many physical systems in which high densities are associated to high temperature. The 



quantification of the amount of order exhibited by the system is done by means of a so- 
called order parameter, which usually ranges between and 1 and which increases with 
the amount of order in the system. 

When the density or noise intensity are varied, the order parameter exhibits a behavior 
which is similar to those of phase transitions in physics (see e.g. [25| [26]). Two typical 
behaviors can be observed. The first one is that of first-order (or discontinuous) phase 
transitions: in this case, there is a parameter range of metastability in which the disordered 
state and the ordered one coexist. The transition from disorder to order or vice- versa 
results in a jump of the order parameter. Additionally, the jumps in either ways do not 
occur for the same value of the density (or noise), leading to a hysteresis behavior. The 
second behavior is that of second-order (or continuous) phase transitions. In this case, the 
transition from disorder to order results in a continuous (but singular in its derivative) 
change of the order parameter. In the region of abrupt change, the system is said to be 
in a critical state. 

Critical phenomena associated to self-organizing systems are grouped into the category 
of 'Self-Organized Criticality' |1]. By contrast to physical systems where critical states 
are reached for very particular combination of the thermodynamic parameters, the critical 
states of self-organizing systems appear like attractors of the dynamics, leading to the 
belief that most of living or social systems operate at the critical state. 

We refer the reader to [57] for a review of these subjects. The modeling of such self- 
organizing systems offers a number of new mathematical challenges, some of which are 
reviewed in the present work. One of the most important challenges, as discussed in 
[57] is the lack of conservation relations, which are the corner-stones of the macroscopic 
theory of large particle systems in physics. Another challenge is related to the possible 
breakdown of the propagation of chaos property, which is another corner stone of the 
statistical mechanics theory of large particle systems. 

As discussed above, these systems exhibit phase transitions. We have emphasized 
the importance of symmetry-breaking phase transitions, which bring the system from 
disordered to ordered states. The broken symmetry is most of the time the rotational in- 
variance symmetry [25112^ . as collective motion exhibits coherent motion in one preferred 
direction. Most often, this preferred direction is random and emerges spontaneously from 
the breaking of the rotational invariance of the fully disordered or isotropic state. But 
another kind of phase transition occurs, and is related to the finite size of the particles. 
It is the transition from the unjammed state (where the particles have room to move 
independently one from each other) to the jammed state, where they are at contact to 
each other, and any motion of one of them induces correlated motion of the neighbor- 
ing particles, sometimes over very large distances. This jamming transition also calls for 
specific mathematical techniques [35] . 

In this work, we will discuss the first challenge, i.e. how to obtain macroscopic models 
for systems which lack conservation relations. We will also briefly review questions related 
to symmetry-breaking phase transitions (a more detailed review can be found in [25]). We 
will not discuss the breakdown of propagation of chaos, and refer the interested reader to 
[T3| [T3|. We will also leave the discussion of the jamming transition for a future review. 
We just mention that this question has been theoretically discussed in relation to traffic 
models in [2], El [HI [23] and to herding in [33]. It has been numerically investigated in 



[271 Eg. 

This review is based on a sequence of papers about the derivation of hydrodynamic 
models for non momentum-conservative particle systems [301 [321 [3S1 HH] and to a lesser 
extent, on [3 [HI [23 I2S1 123 1331 1311 I3S1 133 • There is a vast literature on the mathematical 
modeling of collective motion and self-organization. We refer the reader, e.g. to [H [6l [101 
[m [191 [17] and to the review |57]. 

The paper is organized as follows. In section [2l we introduce the self-propelled particle 
dynamical system which will be at the heart of the present work. Then, in the limit of 
a large number of interacting particles, a mean-field model can be introduced and scaled 
in section [31 The hydrodynamic limit of the scaled mean-field model is studied in section 
[H and gives rise to the Self- Organized Hydrodynamic (SOH) model. Some properties 
of the SOH model are described in section [51 Local existence of smooth solution and 
characterization of weak solutions are investigated in section [SI A discussion of the model 
is given in section [3 and a conclusion is drawn in section [HI 

2 Self-propelled particles interacting through align- 
ment 



Our starting point is the Vicsek model [56] . It describes self-propelled particles modelled 
as particles with constant velocity. They interact with each other by aligning with their 
neighbours up to a certain noise. Originally, the Vicsek model is a time discrete model, 
defined at discrete times separated by a constant time interval At. The positions and 
velocities of A^ individuals at time t" = nAt are denoted by {X]^)k=i,...,N and {y,^)k=i,...,N- 
We denote by m the spatial dimension, i.e. XJ! G M"* (in practice, m = 2 or 3). We 
assume that the particle speeds are constant and uniform, equal to c > 0. Therefore, the 
velocity can be written V^ = cV^, where V^" lies on the sphere S"^~^. The positions and 
velocities are evolved according to the following discrete algorithm: 

XJ:+' = XJ: + cV,- At, (2.1) 

Vfc" = TZ^^^^enV,^, (2.2) 

Eq. (12. ip defines how the position of particle k is updated from time t^ to time t"+^. It 
consists of a simple Euler discretization of the relation defining the velocity as the time- 
derivative of the position. Eq. (12. 2 p states that the new particle velocity is set to the 
average neighbors' direction V^ up to a random noise, expressed by the operator T^^n gn. 
As expressed in (12.31) . V^" is given by normalizing the vector obtained as the sum of the 
velocities VJ^ of the particles lying in a ball or radius R around the subject's position XJ! 
(see Fig. [T]). The quantity R is the interaction range of the subjects. For 6' G M and 
w G §^~^ where S^~^ is the (m — 2)-dimensional unit sphere of the hyperplane orthogonal 
to V, the operator 7lw,eV is defined by 

n,,,0V = cose V + sine w. (2.4) 



In ( I2.2p . w^ are independent, uniformly distributed random vectors on the sphere S^„ ^ 

k 

and 6]^ are independent, uniformly distributed random numbers in some interval [0,D], 
with 0< D <7r. 




Figure 1: The neighbors of the k-th particle located at Xk with velocity Vk (red arrow) 
are found in the ball enclosed by the red circle (ball centered at X^ with radius R) and 
have velocities indicated by blue arrows. The average direction to which particle k aligns 
(up to noise) is that of the sum of the blue vectors. 



In this model, the time step At plays two different roles. The first one is that of a 
time discretization parameter (see Eq. (12. ip ). The second one is that of an interaction 
frequency. Indeed, particles align with their neighbors at each time step. If At is reduced 
by a factor 2, the particles interact twice more frequently. Therefore, the limit At — )■ 
in (12. ip . (12. 3p does not lead to a well-posed problem because the number of interaction 
becomes infinite in a finite time interval. In order to bypass this problem, we decouple the 
two time scales and we introduce a collision frequency u which is independent of the time 
step At. The original Vicsek model is recovered when i/ = 1/At. The modified discrete 
problem when u ^ 1/At is explicitly written in [32] and is omitted here. It can also be 
found in [22] . We directly introduce the time continuous Vicsek model, which is obtained 
as the limit At — )■ of the modified time-discrete problem with collision frequency u. 

Let Xk(t) G M™ and Vfc(t) G S™~^ be the position and velocity of the k-th particle at 
time t. The time-continuous version of the Vicsek model is written as follows. 

Mt) = cVk{t), (2.5) 

dVk{t) = Pyi^ o {u Vkdt + V2DdB^), (2.6) 

Vk = ^^, Jk= Y. ^^■' (2.7) 

where, for V G S™"^, Py^ = Id — V^ (g) ^ is the orthogonal projection onto the plane 
orthogonal to V. Eq. (12. 5p is the formal limit At — )■ of (12. ip . while eq. (12. 7p is the 
same as (12. 3p . The main change is Eq. (12. 6p for the evolution of Vk- It takes the form 
of a Stochastic Differential Equation (SDE). The projection operator Py± ensures that 
the resulting solution Vk{t) stays on the unit sphere, provided that the SDE is taken in 



the Stratonovitch sense (which is indicated by the symbol o). The first term inside the 
bracket is the interaction term. It corresponds to a force acting in the direction of Vk of 
intensity u. The second term is a white noise consisting of independent Brownians B^ in 

113)2 „f ;vif^vic,i + ,r , /O n Tln^ for^f flnof fln^ T^vr^ i ^^f i r^vi P . ^ ^ R^ 



of intensity y2D. The fact that the projection Pyx o dB^ gives rise to a Brownian 
motion on the sphere provided that the SDE is taken in the Stratonovich sense can be 
found in e.g. [H]. We can recover the original Vicsek model through a time discretization 
such that vl\t = 1 and the replacement of the white noise by the process described at 
(12. 4p [22] • The construction of the force term dVk{t) is illustrated in Fig. [2] 



v,dt y^dBt 




Figure 2: Construction of the force term dVk{t) in dimension m = 2: The velocity Vk{t) 
(blue arrow) and the average neighbors' velocity Vk{t) (red arrow) are both vectors of 
the unit sphere §^ (black circle). To define the increment dVk{t), we add the interaction 
force vector uVkit) dt (green arrow) and a small random vector in M^ (yellow arrow) and 
project the resulting vector onto the line normal to Vk{t) (purple arrow). 



Now, letting the number of particles A^ — !■ oo, a mean-field model is obtained. This 
model is described in the following section. 

3 Mean-field model and scaling 

The mean-field model describes the evolution of the one-particle distribution function 
f{x,v,t), which depends on position x G M"*, velocity v G S*""^ and time t > 0. The 
model is written as follows: 

dtf + cv ■ VJ = -V. • (Fff) + DAJ, (3.1) 

Ffix,t) = uP,.Vfix,t), ^^(x,t) = ^M_, (3.2) 



Jf{x,t)= K{ — - — ) f{y,w,t)wdwdy, (3.3) 

Eq. (13. ip is a Fokker-Planck equation. The left-hand side expresses the rate of change of 
/ due to the spatial transport of the particles with velocity cv while the first term at the 
right-hand side takes care of the transport in velocity space due to the interaction force 
Ff. Finally, the last term at the right-hand side is a velocity diffusion term which arises 
as a consequence of the Brownian noise in particle velocities. Because v lies on the sphere 
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gm-i^ Aj,/ stands for the Laplace-Beltrami operator on the sphere. The velocity diffusion 
coefficient D is related to the Brownian noise intensity V2D acting on the particles. The 
force term is proportional to the average neighbors' direction Vf{x,t) around x at time t, 
projected on the hyperplane normal to v through the action of P^±. The proportionality 
coefficient is nothing but the interaction frequency u. The projection P^± ensures that 
the force term is normal to v, as it should to be consistent with the fact that v G S*""^. 
The symbol V^- stands for the divergence of tangent vector fields to §™^^ (and later 
on Vjj will denote gradients of scalar fields on S"^"-*^; we recall that A^ = V^ ■ V^). 

The local average neighbors' direction Vf{x,t) is obtained through the normalization 
of the vector J'f{x,t). This vector is computed by averaging the neighbors' velocities w 
weighted by a function K(^ '^~^^' ^ depending on the distance between the particle and its 
neighbor |x — y|. The average is taken over the probability density f{y, w, t) dydw. In the 
time-discrete case (section [2]), we always took an abrupt cut-off of the interaction region, 
meaning that K{C,) is the indicator function of the interval [0,1]. However, smoother 
cut-offs represented by generic functions K can be taken. The key assumptions on K is 
that it should be positive, integrable, with a finite second moment. The quantity R > 
describes the typical interaction range. The derivation of the mean-field model 03.11] - (13.31) 
from the discrete system (12. 5p - (12. 71) has been performed in a slightly different framework 

inp. 

Macroscopic models are intended to describe the system at large time and space scales 
compared to those attached to the individuals and their interactions. In order to highlight 
the role of the time and space scales, we first non-dimensionalize the mean-field model 
(13. ip - (13. 3p . We let to be a time unit and xq = cto, /o = 1/x™, Fq = 1/to be the associated 
space, distribution function and force units. We introduce the scaled variables x'; t', etc. 
by the following relations: 

x = xox', t = tot', f = fof', F = FoF'. 

Changing from variables (x, v) to {x', v') in system (13. ip - (13. 3p leads to the following model 
(dropping the primes for simplicity): 

dtf + v- VJ = -V. ■ (Fff) + DA J, (3.4) 

Ff{x,t) = iyP^±Vf{x,t), Vf{x,t) = -^—^ (3.5) 

P I I 

Jf{x,t) = K[ ) f{y,w,t)wdwdy, (3.6) 

J{j/,to)GR'"x§™-i V 

where z/ = uto, D = Dto, r] = R/xq. 

Now, we assume that at the chosen time scale to, we have u = O(-) and D = 0(-), 
where e <^ 1 is a small parameter. Specifically, we let: 

z> = -, - = d = 0(l). (3.7) 

e V 

The parameter e measures the interaction time and interaction mean free path, i.e. the 
time and distance needed by a particle to make a finite change in direction of motion due 
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to the interaction force. Equivalently, because the interaction force and noise operators 
are of the same order of magnitude thanks to (13.71) . e is also the time needed by a particle 
to make a finite change of direction due to velocity diffusion. The time and space units to 
and Xq are macroscopic ones while the interaction time and mean free path are microscopic 
quantities. During a macroscopic time unit to? there are a large (i.e. 0{l/e)) number 
of interactions and diffusions that sum up and contribute to making the corresponding 
operators large (specifically, both are 0{l/e) compared to the left-hand side of (I3.4p ). 
At the macroscopic scale, the interaction force and diffusion terms must almost cancel 
each other to yield the order 1 term at the left-hand side of (13.41) . This cancellation is 
the mechanism which forms the 'Local Thermodynamical Equilibrium' detailed below in 
section HI 

We also assume that ?7 ^ 1 and we consider two scaling laws for tj: 

(3.8) 
with r/o > 0, ?7o = 0{l). (3.9) 

In the first scaling (13. 8p . the interaction range is of the same order as the interaction 
mean free path and is therefore microscopic. With the second assumption (13. 9p . the 
interaction range rj = 0{^/e) is large compared to the interaction mean free path. It 
means that a particle interact with many more particles than just those that they are 
able to reach within an interaction time. In the macroscopic limit e — ?• 0, the interaction 
range shrinks to 0, meaning that it is smaller than a macroscopic quantity. Therefore, 
with the second scaling, the interaction range is intermediate between the microscopic 
scale and the macroscopic one. This choice of scales is illustrated in Fig. [31 



MacrS\scale = 0(1) 

IntcractionVangc = 0{\/s) 

• •* * 1 

tion mcanlfrcc path = 0{e) 





Figure 3: The different scales of the problem in the case (ii) (see (13.91) ): the microscopic 
scale is the interaction mean free path and is 0{e). The interaction range is the in- 
termediate scale and is 0{y/e). Finally the macroscopic scale is the scale of the whole 
system. 



With these scahng assumptions, the scaled mean-field model is written: 

e{dtr + V ■ V.r) = -V„ ■ {F]J') + dAJ', (3.10) 

F]{x,t) = P,,.v]{x,t), v]{x,t) = ^f^, (3.11) 

\Jf [x,t)\ 

J]{x,t)= I Ki}y—^)f{y,w,t)wdwdy, (3.12) 

-'(i/,u;)gR'"x§'"-i 1 

The goal is to study the limit £ — ;■ 0. Since 77 — ;■ in both scalings (13. 8p and (13. 9p . we can 
expand F^{x, t) in powers of 77 and get [30] : 

F]{x,t) = F;{x,t)+r^'F}{x,t) + 0{v'), (3.13) 

F;{x,t) = P,^Ufix,t), (3.14) 

F/(a:,t) = Ap^xP„xA,,j^, (3.15) 

where the local density pf, local current density jf, and local average direction uj are 
defined by: 



pf{x,t) = / f{y,w,t)dw, 



w& 



jf{x,t)= f{y,w,t)wdw, (3.16) 

Uf{x,t) = ^Ip^ (3.17) 



We have denoted by k the second moment of K, i.e. 

1 



2„, mnKi'di, 



provided that K is normalized to 1 (without loss of generality), i.e. /^gj^m -^'^(1^1) d^ = 1. 
Neglecting terms of order 77^, and defining the interaction (or collision) operator Q{f) 
by: 

Qif) = -V. ■ (P^^Uff) + dAJ, (3.18) 

we can write the system 

e{dtr + V ■ V.r) + r/^V. ■ {Fjj') = Qif), (3.19) 

Now, there are two cases according to which scaling assumption hypothesis (13. Sp or (13. 9p 
is made. In case (13. 8|) . the last term of the left-hand side of (I3.19P is of order e"^ and can 
be neglected. In case (13. 9p . we have 77^ = r/o£ and the last term of the left-hand side is of 
the same order as the other terms of the left-hand side. We can collect them and get the 
problem 

e{dtr + V ■ v.r + ^ov. ■ (F/.r)) = Q(r), (3.20) 

We note that we recover case (13.81) by just letting t/q = in (I3.20p . Therefore, case (13.81) 
becomes a sub-case of case (13.91) and we now only consider (13.201) below. In the next 
section, we give the formal e — )■ limit of this model. 
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4 Self-Organized Hydrodynamics (SOH) 

The formal macroscopic limit e — )■ of ( I3.20p has been studied in [301 132] • We first state 
what are the equilibria, i.e. the solutions of Q{f) = 0. These solutions are important 
because /° = lim£_>.o f^ (if it exists) is obviously an equilibrium thanks to f l3.20p . Since Q 
only operates on the velocity variable, we first focus on the velocity dependence of these 
equilibria. Denote by £" = {/ 1 Q{f) = 0} the set of equilibria. Then, it is possible to 
show that, under reasonable regularity assumptions, S is spanned by the so-called von 
Mises-Fisher (VMF) distributions, i.e. 

S = {v -^ pMu{v) for arbitrary peR+, ue S"""^}, (4.1) 

where 

Muiv) = ZJ^ exp{——), Zd= exp{——)dv. 



We note that Mu{v) is a probability density and that Z^ does not depend on u. The 
element u G S"^"^ is called the direction o 
concentration parameter. We also note that 



element u G S"^ ^ is called the direction of the VMF distribution while k = ^ is its 



v& 



f ^sm-1 exp (^) (u ■ v) dv 
MJv)vdv = ciu, ci = ci(d) = -^"^^r ^ '';\ • 4.2 



The quantity ci{d) does not depend on u and satisfies < ci{d) < 1. Eq. (14. 2 p shows 
that the current associated to M„ is directed and oriented along u and its magnitude 
is defined by Ci. The function k = 1/d G [0, oo) i— )■ Ci{d) G [0,1] is strictly increasing, 
onto. Small values of Ci correspond to VMF distributions close to the uniform isotropic 
distribution, while values of ci close to 1 correspond to VMF distributions close to Dirac 
deltas. Therefore, the parameter ci measures the degree of alignment of an ensemble 
of particles whose orientation v is statistically defined by M^. It is used as an order 
parameter in the study of phase transitions between disordered and aligned states. Fig. |4] 
depicts a VMF distribution in dimension 1. 

The VMF distribution allows us to write the collision operator Q in the form of a 
Fokker-Planck operator. We have: 



Q{f){v) = dV,- 
We deduce a dissipation estimate 



Muf{v)V, ' 



M„, 



/ Qif)j^dv = -d[ 



f ' 



M„^. dv <0, 



with an equality if and only if / is an equilibrium i.e. f ^ S. The proof of (14. ip relies 
on this inequality. Unfortunately, this dissipation estimate does not yield a H-theorem in 
the spatially inhomogeneous case because the multiplier M„ depends on moments of /. 
From the determination of the equilibria, we know that /° = lime_j,o /^ (if it exists) is 
such that 

f{x,v,t) = p{x,t)Mui.,t){v), (4.3) 
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Figure 4: The VMF distribution in dimension m — 1 = 1, i.e. when the velocity variable is 
defined as w = e*^ with x g] — tt, it] (illustration taken from Wikipedia). The concentration 
parameter 1^ = 2 controls the width of the distribution. For k = 0, i.e. d — )■ 00, the 
VMF distribution is very close to a uniform distribution (black horizontal line). When k 
increases or d decreases, the height of the bump increases and its width decreases: the 
blue, blue-green, green, red and purple curves correspond to increasing values of k or 
decreasing values of d. The direction of the VMF corresponds to the center of the peak. 
Here, the direction of the VMF is x = 0. 



where now p = p{x,t) and u = u{x,t) may have non-trivial dependences upon (x,t). 
Indeed, since Q operates only on the v variable, it does not impose any constraint on the 
dependence of the parameters p and u of the equilibrium on {x, t). 

In order to find this dependence, we need to introduce the second important concept, 
which is that of a collision invariant. A classical collision invariant is a function ip{v) such 
that 

/ Q{f){v)tlj{v)dv = 0, 

for all distribution functions / with reasonable regularity (we will not dwell on this point 
here). The set of these collision invariants is denoted by C and is a vector space. The 
parameters (p, u) G M+ x S'^^^ of the equilibrium span a space of dimension m. To 
specify them, we need m independent collision invariants, i.e. C must be of dimension m. 
Collision invariants are strongly related to conservation laws. For instance, in classical 
gas dynamics, the conservations of mass, momentum and energy provide m + 2 indepen- 
dent collision invariants which allow to determine the m + 2 independent parameters of 
the equilibrium Maxwellian distribution. In the present case, there exists only one con- 
servation relation, namely the conservation of mass, which is associated to the collision 
invariant ip = 1 (indeed, by Stokes formula, Eq. 03.18P immediately leads to the fact that 
J^gg^.i Q{f) dv = 0). Due to the self-propelled character of the particles, there exists no 
momentum nor energy conservation nor any other kind of conservation. Consequently 
there exists no collision invariant independent from ip{v) = 1. The dimension of C is 1 
and is strictly less that m. So, m — 1 independent collision invariants are lacking to fully 
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determine the parameters p and u of the VMF equihbria. 

To overcome this problem, we need to weaken the concept of a colhsion invariant. 
We define the concept of a 'Generahzed Colhsion Invariant' (GCI) as follows. For any 
arbitrary vector u G S™~^, we first define the linear operator Q„ by 



QuU){v) = dV, 






we obviously have Q{f) = Qufif)- Then, for any m G S™ ^, a function ip^ is a GCI 
associated to u if and only if we have: 

/ Qu{f){v)Mv)dv = 0, V/ such that P„xM/ = 0. (4.4) 

In other words, instead of testing ip against Q{f) for all /, we restrict the set of test 
functions / to those whose mean velocity direction Uf is proportional to u (and since 
both are unit vectors, this means Uf = ±n) and use the linear operator Q^ instead of the 
nonlinear one Q. By taking a smaller set of test functions /, we are enlarging the set of 
possible solutions ipu- The concept of a GCI is now associated to a choice of a direction 
u. By a duality argument based on the fact that both the operator Qu and the constraint 
Pu-^Uf = (which is equivalent to Pu^jf = 0) are linear in /, it can be shown [32] that 
04.4p is equivalent to the existence of a vector f3 G M™ with (3 ■ u = such that 

Ql{^^){v)=(3-v, (4.5) 

where Q* is the L^-adjoint of Q„. It is easy to see that the problem of finding ipu is linear 
and consequently, the set Qu of the GCI associated to u is a linear vector space. Since Q„ 
is a second order elliptic operator, Eq. 04. 5 1) is a second order elliptic problem posed on 
S'""! whose solution can be found using Lax-Milgram theorem. This leads [211 [261 EH] to 
the fact that Qu is an ?ti- dimensional vector space given by 

g^ = {y ^ C + h{u ■v)(5- V, with arbitrary C G M and /3 G M'" with /3 • n = 0}. 

The unique scalar function h is defined as follows. Denote hj u ■ v = cos 6. Then, 
h{cos9) = g{9)/ sva.9 where g{6) is the unique solution in the space V of the elliptic 
problem L*g = sm6, with 

and 

V = {g \ (m-2)(sin^)f"2^GL2(0,7r), (sin^)f "^ ^ G i/o'(0, vr)}, 

(denoting by Hq{0, it) the Sobolev space of functions which are square integrable as well as 
their derivative and which vanish at the boundary). Since a function of Qu is determined 
by the scalar C and the vector /3 such that /3 ■ u = 0, it depends on m independent 
parameters and forms a linear space of dimension m. 

Now, we can perform the formal limit e -> 0, assuming that /^ — > /° as smoothly 
as needed (i.e. with a regularity such that all the formal computations below can be 
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rigorously justified). We already know that f^ is an equilibrium given by (14.31) . Now, 
integrating (13.201) with respect to v, using Stokes' formula on the sphere and taking the 
limit £ —7- 0, we find the continuity equation 

dtp + Vx ■ (cipw) = 0. 

The fact that the particle flux at equilibrium is expressed by cipu comes from ( 14. 2p . 

In order to find the equation satisfied by u, we multiply (13.201) by the vector-valued 
GCI "^Ufe, where, for an arbitrary vector u G S^'^, \l/«(f) = h{u ■ v)Pu±v, and take 
the limit e — ?■ 0. We note that, as a consequence of the definition of the GCI, we have 
/^g§m-i Qif^) ^«/e dv = 0. Consequently, the right-hand side of (I3.20p . which is singular 
as 5 — )■ 0, disappears after integration against "^Ufs- After some algebra, the left-hand 
side gives rise to the following equation for u: 

p {dtu + C2(m ■ Vx)u) + dPu±VxP = C3Pu±A^{pu) 

with 

C cos 9 h{cos 9) e^ sin"" 9 d9 
Jq h{cos 9)e d sm 9 d9 

We summarize the result in the following: 

Theorem 4.1 As e ^ 0, we have f^ — ?■ /° formally, where f^ is an equilibrium ( [^.3p 
whose parameters p{x,t) andu{x,t) satisfy the following 'Self-Organized Hydrodynamics' 
(SOH) system: 

dtp + V, ■ (cipn) = 0, (4.6) 

p {dtU + C2{u ■ V^.)w) + dP^±V^p = C3Pu^A^{pu), (4.7) 

\u\ = 1. (4.8) 



In the next section, we review some properties of the SOH model, we state the available 
existence results and we address the question of its numerical approximation. 

5 Properties of the SOH model 

The SOH system (14. 6 1) - (14. 8 p bears strong similarities with the isothermal compressible 
Navier-Stokes (NS) model, as recalled below: 

dtp + V.-{pu) = 0, (5.1) 

p {dtu + {u ■ V.^)u) + TV^p = uV^- (V^u + V^u*), (5.2) 

where p{x,t) > and u{x,t) G M'" are the gas density and mean velocity, T is the 
temperature (the analog of our d) and u is the viscosity (the analog of our C3). The 
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quantity VxU is the tensor gradient of the vector field u and the exponent '*' denotes the 
transpose of a tensor. We adopt the notation V^- indifferently for divergence of vectors 
and tensors. We have written the momentum equation fl5.2p in non- conservative form to 
highlight the analogy with fl4.7p but it is more natural to express it in conservative form. 
This form is obtained by multiplying (15. ip by u and adding to (15.21) . This leads to 

dt{pu) + V^ ■ {pu ®u)+ Vx{pT) = Vx ■ (z/( V^u + V^u* )) , (5.3) 

where u® u stands for the tensor product of the vector u by itself. In the inviscid case 
1/ = 0, the conservative form allows to define solutions with jumps thanks to the Rankine- 
Hugoniot condition. Additionally, the entropy condition, which we will not recall here, 
permits the selection of a single solution among the possible weak solutions. 

Both the SOH and NS models are constructed as nonlinear systems of first order 
equations perturbed by diffusion. The nonlinear first order parts consist of the left- 
hand sides of (14. 6p . (14. 7p on the one-hand and (15. ip . (15. 2p on the other hand. The 
diffusion operators appear in the velocity equations (the right-hand sides of (14. 7p and 
(15. 2 p respectively). Additionally, both nonlinear first order parts are hyperbolic. This 
is standard for the NS model: its first order part is the isothermal compressible Euler 
equations whose characteristic velocities are A± = (m ■ (^) ± vT and \q = u ■ ^ where ^ is 
the unit vector in the direction of propagation of the characteristic wave. In the case of 
the SOH model, these characteristic velocities have been computed in |32] and are given 
by 
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1 r, . . . /, .o , ..o ..o. \ 1/2 



(ci + C2)u-i± ((C2 - cif {u ■ if + 4rf (1 - (n ■ if)) ' , Ao = ci n ■ e (5.4) 



The first important difference between the two models is that the SOH model involves 
the geometric constraint (14. 8 p which imposes the velocity u to be of unit norm. By 
contrast, the standard NS model does not involve such a constraint. We note that, for 
smooth solutions of the SOH model, this geometric constraint is satisfied at all times 
provided it is satisfied initially, i.e. if | u\t=Q | = 1- Indeed, P„± multiplies the last term of 
the left-hand side and the right-hand side of (14.71) . Therefore, multiplying scalarly (14. 7p 
by u, we get: 

hyp 
[dt + C2(n ■ VJ) i^-j-) = [{dt + C2(n ■ V.,))u] ■ n = 0, 

which shows that \u\^ is a conserved quantity along the characteristics of the vector field 
C2u{x,t). If |mP is initially equal to 1 uniformly, it stays equal to 1 uniformly at future 
times. This computation requires the solution to be smooth. The preservation of the 
norm along discontinuous trajectories is still a conjecture at this time. 

The geometric structure is brought to the model by the multiplication of the last term 
of the left-hand side and of the right-hand side of (14.71) by the projection operator P^± . The 
matrix P^± = Id — u ®u being a non-trivial function of u, this multiplication introduces 
non-conservative products. Indeed, the corresponding terms are spatial derivatives of 
non-trivial functions of the unknowns p and u, multiplied by non-trivial functions of these 
unknowns. Non- conservative hyperbolic systems are reviewed in [16]. They have some 
unpleasant features such as the impossibility of defining shock relations for discontinuous 
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solutions by means of the Rankine-Hugoniot condition. As a consequence, the shock 
speeds are not well-defined and there is no simple criterion to single out a particular 
solution among the possible ones. 

In the SOH model, the fiuid mean velocity is not given by u but rather, by Ciu. The 
vector u gives the average orientation of the particles. Most continuum models of self- 
propelled particles available in the literature, such as [531 EH ES] , use the mean velocity 
as a variable of the model, instead of the mean orientation. As a consequence, in these 
models, u is not constrained to be of unit norm, and they bear a closer resemblance with 
the NS than the SOH model presented here. However, these models are constructed on 
phenomenological bases, while the SOH model is derived through a rigorous hydrodynamic 
limit of the underlying mean-field model. 

The second important difference between the SOH model and the standard NS model 
is the constant C2 which, in general, is different from Ci. In the present case, we can prove 
[32] that < C2 < Ci (where C2 = Ci if and only if the noise intensity d = 0). Since CiU is the 
flow velocity and C2U is the velocity at which the velocity itself is transported, this means 
that velocity is transported upstream the flow. The transport of velocity is a transport 
of the information about how agents should adjust their velocity to accommodate for the 
presence of other agents in the front. This is a situation similar to vehicular traffic (see 
e.g. [3]), where information about how drivers should adapt their velocity to the vehicles 
in the front propagates upstream the fiow and even, in congestion situations, opposite to 
the fiow direction. 

When anisotropic vision is considered, like in [38], arbitrary values of C2 G M can be 
generated. Anisotropic vision means that the interaction kernel K does not only depend 
on the distance \y — x\ between the subject located at x and its partner located an y 

but also on the bearing angle, i.e. the angle {v,y — x) between the subject's velocity v 
and the line of sight y — x under which he sees his partner. In particular, forward vision 
corresponds to situations where K is lower for obtuse bearing angles, and backwards 
vision, where K is larger. In the isotropic vision case discussed so far, K only depended 
on the distance \y — x\. In [38], it is shown that forward vision results in general in 
lower values of C2 than isotropic vision, while backwards vision results in larger values. In 
the case of forward vision, for large enough values of the noise intensity d, one can even 
have C2 < 0, meaning that velocity is transported in opposite direction to the fiow. As 
noted above, this is similar to car traffic in congested situations. In the case of backwards 
vision, for large values of d, we can have C2 > Ci, meaning that information propagates 
downstream the fiow. It is shown in [27] that C2 < Ci increases the ability of the model 
to generate density concentrations and shocks, while this tendency is lowered in the case 
C2 > ci. As pointed out above, the case C2 < ci bears analogies with vehicular traffic. 
By contrast, the case C2 > ci is exemplified in locusts' coordinated mass migration. The 
locust species Schistocerca gregaria is a cannibalistic species and individuals try to avoid 
being bitten by others in their back [7]. This is an example where backwards vision (in 
this case, the sensory mechanism is not vision but rather abdomen innervation) controls 
collective motion. 

The SOH model is not Galilean-invariant because of the constraint |m| = 1. Indeed, it is 
not possible to translate the velocity m by a uniform velocity V and keep this constraint. 
That C2 7^ Ci also induces some non-Galilean effects. Indeed, if we perform the same 
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manipulations as those leading to the conservative form (I5.3P of the NS model, we are led 
to: 

dt{pu) + V^ ■ {cipu (g) n) + (c2 - ci) {u ■ V^)u + P^±Vr,{dp) = cgP^x A^.(pn), (5.5) 

The extra term (c2 — ci) {u-'Vx)u can be understood as a velocity-dependent extra pressure 
term. The fact that pressure is velocity dependent is a signature of a non-Galilean model, 
because in Galilean-invariant fluids, the pressure is an intrinsic property of the fluid, which 
is the same in all reference frames. Another manifestation of the non Galilean-invariant 
character of the model is the expression of the sound speed Cg, which corresponds to the 
square root term in the expression (15. 4p of the characteristic speeds A±, i.e. 

c,= {{c2-c,f{u-0' + ^d{l-{u-0')Y^'- (5.6) 

This expression depends on the angle between the propagation direction ^ and the ve- 
locity u. In Galilean-invariant fluids, the sound speed is independent of the propagation 
direction. Like our model, the previously proposed continuum models of swarming are 
not Galilean-invariant and the effect of non Galilean-invariance on sound propagation is 
reviewed in [55] . 

Finally, we note that special stationary solutions are given by mills, i.e. solutions of 
the form 

p{x) = pi\x\) = po{r /roY^^'^, u = x^/\x\, 

where po > and tq > are arbitrary constants. In such solutions, the average velocity 
corresponds to a spinning about the origin at a constant speed, while the density profile is 
unbounded and increases like a power law of the distance to the origin. According to the 
position of C2/d with respect to 1, the density profile is convex or concave. Convex density 
profiles are associated to small noise and consequently, to a clear spinning organization of 
the individual particle velocities. By contrast, at large noise, the density profile becomes 
concave. Due to the large noise, the individual velocities show large deviations around 
the average spinning velocity, and the organization of the individual particles in a milling 
pattern is less clear. These are illustrated in Fig. O The question of the stability of these 
milling solutions is an open problem. Milling is a common yet intriguing social behavior 
in fish populations 



6 Local existence of smooth solutions and selection 
criterion for discontinuous solutions 



In this section we state the existence results of [30]. We first start with the case of 
dimension m = 2. In this case, we can take c^ > 0. We assume that the spatial domain is 
the two dimensional torus, i.e. x G 11^, with 11 = [0, 1] with periodic boundary conditions. 
The initial data {pq.uq) are such that po > and |mo| = 1. We write u = (cos (p, sin (p). 
We have the 

Theorem 6.1 IWf We assume that the initial data belong to H^(lP) with s > 2. Then, 
there exists a time T > and a unique solution {p,ip) G L°°{[0,T],H^{lP)) fl H^{[0,T], 
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Figure 5: Schematics of a milling solution: left: mean orientation u (blue arrows) and 
individual particle velocities Vj!^ (green arrows) ; right: density profile p as a function of 
the distance to the origin r = \x\. Top figure: small noise case C2/d > 1: the density 
is convex and the particle velocities V^ follow closely the mean orientation u. Bottom 
figure: large noise case C2/d < 1: the density is concave and the particle velocities V^ 
deviate strongly from the mean orientation u. 
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if" ^(n^)) of the SOH model such that p remains positive. If, in addition, C3 > 0, then, 
the solution also belongs to L'^{[0,T], H'+^iW^)). 

The proof relies on the fact that the first order part of the SOH model admits a 
symmetrizer. This symmetrizer allows the development of energy estimates. 

We now turn to the dimension m = 3. In this case, we use global spherical coordinates 
and write u = [sin 6 cosy?, sin 6' simp, cos 6). Due to the singularity of this coordinate 
system for 6^ = and 6^ = vr, we choose initial velocities such that (6*0,930) ^ [(^m,(^M] x 
[0,27r], with < 9m < Om < ^i- In the case C3 = 0, we can use the finite speed of 
propagation of hyperbolic models to guarantee that the solution 6 will stay away from the 
singular points 6 = and 6 = tt during some interval of time. Because diffusion involves 
propagation at infinite speed, this property is lost in the case C3 > and the proof 
breaks down. This is why we restrict ourselves to the inviscid case C3 = in dimension 
m = 3. Most likely, this restriction could be waived at the expense of technicalities, 
such as changing the local coordinates on the sphere near the singular points but these 
developments are left to future work. 

Again, a periodic domain 11^ is chosen as spatial domain and the initial density is 
supposed positive po > 0. Then, in the case m = 3 and C3 = 0, the theorem reads as 
follows 

Theorem 6.2 fSUf We assume that the initial data (po,^0)V^o) belong to H^{lP) with 
s > 5/2 and that po > 0, sin^o > 0. Then, there exists a time T > and a unique 
solution {p,e,ip) e L°°{[0,T], H'(n.^)) n H\[0,T], H'-\U^)) of the SOH model such that 
p remains positive. 

We show that the inviscid SOH model in dimension m = 3 can be written in the form 
of symmetrizable hyperbolic system. Existence and uniqueness follows from the classical 
theory of symmetrizable hyperbolic systems. 

These local existence results for smooth solutions do not give any information about 
the existence and uniqueness of non-smooth solutions. As a non- conservative model, the 
SOH may have multiple shock velocities [16] and there is no obvious analytic criterion 
to single out one particular solution. Following [IS], given two states {pi,Ui) and {pr,Ur) 
separated by a spatial discontinuity, there are as many shock relations as possible paths in 
the state space spanned by (p, u) (i.e. IR+ x S*""^) connecting these two states. Most non- 
conservative systems, like the SOH model, are obtained by some coarse graining procedure 
from an underlying 'microscopic' model. What is the correct path is an information which 
has been lost in the coarse graining procedure. In principle, this information could be 
retrieved if one could 'interrogate' the microscopic model. Unfortunately, this is not 
doable in practice. 

Here, the coarse graining procedure is the hydrodynamic limit and the microscopic 
model is the time-continuous particle system (12.51) . (12. 7p . Therefore, we can compare how 
the SOH and particle systems resolve a Riemann problem and try to find experimentally 
what are the correct shock relations for the SOH model. In fact the problem is considerably 
simplified by the remark made in ^\. It is based on the observation that the SOH model 
can be formally obtained as a relaxation limit of a conservative model, the so-called 



Relaxed SOH model (RSOH). The RSOH model is written as follows: 

dtp' + V. ■ {c^p'u') = 0, (6.1) 

dt{p'u^) + V. ■ (C2 p'u^ ® u^) + dV,p' = -\ p' (1 - |M^n u\ (6.2) 



where 5 <^ 1 is a relaxation parameter. Here we have considered the inviscid case C3 = 
for simplicity. It is possible to formally show [19] that (p^, u^) — )■ (p, n) as 5 — ?■ where 
(p, u) is a solution of the SOH model. The RSOH model is hyperbolic if and only if 
C2 > Ci. By contrast, in the case where C2 < Ci, it is only conditionally hyperbolic, when 



d 
u < 






The SOH model being a relaxation limit of the larger conservative system (16. ip . (16. 2p . one 
may think that it is amenable to the relaxation theory for systems of conservation laws 
[T7] . In fact, it is not the case. Indeed, in [17], a key hypothesis is that the relaxed system 
is also a system of conservation laws. Here the system is non-conservative and the theory 
does not apply. In particular, the 0{6) correction terms to the SOH model that occur 
in a Chapman- Enskog expansion of the RSOH system can be computed. They are not 
diffusive terms by contrast to those which are obtained in the theory of [T7]. Indeed, the 
geometric constraint in the SOH model opens a wealth of novel features, most of them 
being still unexplored. 

From the RSOH model, a relaxation scheme for the SOH model can be designed. It 
consists of a time splitting method. Given an approximate solution (p", «"■) of the SOH 
model at time t" = nAt, the first step of the splitting consists in solving the non relaxation 
part of the RSOH model (16. ip . (16. 2p . i.e. the system 

dtp + V, ■ (cipn) = 0, (6.3) 

dt{pu) + \/^,-ic2pu(g)u) + d\/^p = 0, (6.4) 

over one time step At with initial condition (p", m") by a standard shock-capturing scheme 
such as the Rusanov method [5l]. This leads to intermediate values {p"'~^^ , u"'^^) . Then, 
in the second step of the splitting, the relaxation part of the RSOH model is solved over 
one time-step At with initial condition (p"+^,'u""'~^). This relaxation part reads: 

dtp' = 0, 

dt{p'u') = --n'il-\u'\')u'. 


The solution of this system can be explicitly computed and in the limit 6-^0, just reduces 
to the normalization of the velocity issued from the first step of the splitting. This leads 
to the new value (p""'"^,^""'"^) of the solution at time t"+^: 

In the situation where the RSOH model is not hyperbolic, the method can still be 
used. Indeed, in the splitting method, the solution of the non-relaxation part of the 
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RSOH model fl6.3p . (16.41) is immediately reprojected onto a solution of the SOH model 
through the normalization fIG.Sp . So, an instability does not have time to develop before 
the solution is reprojected into a solution of the SOH hyperbolic model. The stability of 
this solution methodology has been experimentally demonstrated in [15] but a rigorous 
proof of this property is still lacking. Now, a pending question still remains: how to prac- 
tically solve the non- relaxation part of the RSOH model (16. 3p . (16. 4p in a non- hyperbolic 
situation ? The answer is easy in the context of e.g. the Rusanov method. Indeed, in this 
method, the computation of the characteristic speeds is only needed in the determination 
of the numerical viscosity. Using the 'wrong' numerical viscosity is not detrimental as 
long as it does not become too small. On the other hand, overestimating the numerical 
viscosity may deteriorate the quality of the solution, but does not threaten the stability 
of the method. In the model, the breakdown of hyperbolicity occurs by the appearance 
of a negative real number inside the square root in the expression of the sound speed 
(15. 6p . Taking the absolute value of this negative number allows to assign a real value 
to the square root and consequently, to set a value to the numerical viscosity. Using 
this value apparently provides correct solutions without the appearance of any instability 
problem [27]. However, a rigorous basis to this methodology is still unavailable. 

It has been experimentally discovered in [12] that using this relaxation scheme to solve 
the Riemann problem for the SOH model provides extremely good approximations of the 
corresponding solution of the particle system. By contrast, any standard shock captur- 
ing method adapted to non-conservative models provides wrong solutions. Therefore, 
although we are not able to provide an analytic criterion for the selection of the 'correct' 
discontinuous solutions, an experimental way to compute them numerically exists. Details 
and numerical evidence can be found in [15] • 

To conclude this section, we present some numerical simulations in Fig. O This 
simulations show that the SOH model provides a fairly good approximation to the solution 
of the particle model even for long simulation times. 
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7 Discussion 

We have presented the derivation of the SOH model in the simplest situation where the 
collision frequency and noise are constant. Several generalizations have been proposed. In 
[32] . the collision frequency v depends on the angle u-Vf between the particle velocity and 
the neighbor's average orientation. All the results given above extend to this case, except 
for the value of the constants Ci, c^ and for the coefficient multiplying the pressure term 
Pu^SJ xP- We denote by 9 this coefficient (in the case discussed here we have Q = d). In 
[32] . the derivation of the SOH model was made in dimension 3 only. In [38], an arbitrary 
dimension m > 2 is considered and both the collision frequency v and the noise d are made 
dependent of the local density p/. Again, all the present results extend to this case, except 
that now, Ci, Ci and 6 become functions of the density. Expressions of these constants 
are given in arbitrary dimensions and asymptotic developments of them when d — )■ and 
rf — ;■ oo are provided. Additionally, it is shown that, with density dependent coefficients, 
the SOH model may lose its hyperbolicity. A detailed analysis of the hyperbolicity of the 
model can be found there. 

In our SOH model the order parameter C\[d) is a constant, since c? is a constant. 
Therefore, it is not able to account for the phase transitions reported e.g. in [SS]. In this 
work, the authors observe that the discrete particle model fl2.1l) - fl2.3p undergoes a phase 
transition from disorder (small values of c\) to alignment (values of c\ close to 1) as either 
the noise intensity d is decreased or the density p is increased. In [T6|, it is observed 
that in some regimes an attractor consisting of travelling bands emerges. In these bands 
the particle density and order parameter are large, indicating a high level of alignment. 
The bands travel through an ocean of low density disordered particles. Therefore, the 
disordered and ordered phase may coexist and are separated by dynamic interfaces. Our 
SOH model with a constant order parameter is unable to account for such features. 

A major step forward to resolving this inaccuracy has been made in a series of works 
[211 [23 [2Sl [3m • In [2H ISH] , the collision frequency v depends linearly on the norm of the 
local current i.e. v = vq\Jj\^ where z/q is a constant. In this case, it is shown that, besides 
the VMF equilibrium, a second equilibrium coexists which consists of the uniform distri- 
bution of orientations. Additionally, at low densities, this is the only stable equilibrium. 
But beyond a critical density pc^ the VMF equilibria emerge as stable equilibria, while 
the uniform equilibria become unstable. Therefore, a second order (or continuous) phase 
transition happens at the critical density p^. In [2S1[2S|, niore general dependences of v 
and d upon |J7/| are considered and a wealth of interesting behavior, such as first order 
(discontinuous) phase transitions with hysteresis are exhibited. Additionally, the analysis 
of [2H [39] is complemented in terms of rigorous rates of convergence of equilibria and 
properties of the macroscopic model. In these works, it is shown that, in the area where 
the density is below the critical density and where consequently the isotropic equilibria are 
stable, the macroscopic dynamics is given by a nonlinear diffusion equation in place of the 
SOH model. The SOH dynamics is recovered above the critical density where the VMF 
distribution exists and is stable. These works pave the way to a rigorous macroscopic 
modeling of the patterns observed in [T6|, [56]. However, to be effective, the theory should 
provide a way to interface the nonlinear diffusion model and the SOH model across the 
interface between the phases. This theory is still missing. 
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Particle model at t=Os SOH model at t=Os 
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Figure 6: Comparison between the particle model f l2.5p . (12. 7p (left panels) and the SOH 
model (right panels). Top panel: time t = 0. Middle panel: time t = 30. Bottom panel: 
time t = 60. The speed of the particles is 1 and the size of the domain is 20. The 
average density in each cell is represented by a color code from low density (light yellow) 
to large density (dark red). The average velocity in each cell is represented by an arrow 
in the center of the cell. For the SOH model, the magnitude of the average velocity is 
constant equal to ci. We observe that the agreement is fairly good in spite of a noisy 
initial condition and a large simulation time. 
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The Vicsek model has been criticized on the basis that the assumption of constant 
speed of motion is unreahstic. It has been argued that in practice, the agent speeds 
in fish schools for instance is not constant. In |3], we consider a model with no such 
constraint on the particle velocities. It consists of the Cucker-Smale (CS) model with 
the addition of a self-propulsion force. The CS model [IS |20l |2ll |l2l Sal |50l |52] has 
received a lot of attention recently. Like the Vicsek model, it describes the relaxation of 
the agents' velocity to that of their neighbors, but there is no constraint on the agents' 
speed. Most previous work on the CS model ignore the self-propulsion velocity. In [S] 
we add the self-propulsion force in the form of an operator whose effect is to relax the 
norm of the velocity to a constant value. This operator takes the form of a Ginzburg- 
Landau type term in the mean-field equation. We perform the successive hydrodynamic 
limit and limit of large self-propulsion force in this model. We show the emergence of 
phase transitions between a nonlinear diffusive regime and a hyperbolic regime described 
by the SOH model. Therefore, the essence of the SOH model lies in the presence of a 
large propulsion force rather than on the normalization of the particle velocities. The 
importance of the self-propulsion force in triggering the phase transition in self-propelled 
particle systems was already remarked in |54] . 

The SOH model also emerges as the hydrodynamic limit of a large class of microscopic 
dynamics. For instance, in [M], it describes the hydrodynamic limit of a system of agents 
controlling their motion by acting on the curvatures of their trajectories (like a driver 
controls the trajectory of his car by acting on the steering wheel) and trying to join their 
neighbor's direction of motion. The Individual model, now referred to as the Persistent 
Turner (PT), has been derived from experimental observations of fish trajectories in |l0]. 
The validation of the model has been made by comparing the diffusion constant induced 
by the model |33] to the experimentally observed one |1D]. In [H], a model close to the 
one proposed in |3l] has been validated by comparisons to experimental data. 

In [29], an extension of the time-continuous particle model has been proposed in di- 
mension 777, = 3. In this model, the interaction force is complemented with a precession 
term. In addition to driving the agent's velocities towards their neighbors' average ve- 
locity, it makes the former rotate about the latter. The resulting hydrodynamic model is 
an extended SOH model with extra transport and diffusion terms acting in the direction 
normal to the average velocity. Interestingly enough, when the self-propulsion speed is 
set to 0, the resulting SOH model is nothing but the Landau-Lifschitz-Gilbert (LLC) 
equations of micromagnetism [12]. Therefore, our theory offers one of the very few (if not 
the only) derivation of the LLC equations from first principles. 

These two examples confirm that the SOH model is a generic model which may apply 
to a large class of situations. Additionally, in [31], the SOH Model provides an example 
of hydrodynamic models which can be derived from the interplay of kinetic theory and 
game theory. 

8 Conclusion 

In this work, we have reviewed the derivation and properties of the Self- Organized Hydro- 
dynamic (SOH) model. It describes the large-scale behavior of systems of self-propelled 
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particles interacting through local alignment. The main difficulty in this derivation is re- 
lated to the lack of conservation properties of the underlying particle dynamics. We have 
described how the use of the Generalized Collision Invariants can overcome this problem. 
The known mathematical properties of the model have been reviewed and its connections 
with other kinds of self-propelled particle systems have been discussed. In the future, 
the spatial dynamics of the interface between the SOH and nonlinear diffusion models in 
the case of phase transitions between ordered and disordered states will be investigated. 
Other kinds of constrained dynamics or geometries will be considered. The question of 
the boundary conditions applying to the SOH model will also be addressed. 



References 

[1] I. Aoki, A simulation study on the schooling mechanism in fish, Bulletin of the Japan 
Society of Scientific Fisheries, 48 (1982) 1081-1088. 

[2] C Appert-Rolland, P. Degond, S. Motsch, Two-way multi-lane traffic model for 
pedestrians in corridors, Netw. Heterog. Media, 6 (2011) 351-381. 

[3] A. Aw, M. Rascle, Resurrection of second order models of traffic flow, SIAM J. Appl. 
Math., 60 (2000) 916-938 

[4] P. Bak, C Tang, K. Wiesenfeld, Self-organized criticality: an explanation of 1// 
noise, Phys. Rev. Lett., 59 (1987) 381-384. 

[5] A. Barbaro, P. Degond, Phase transition and diffusion among socially interacting 
self-propelled agents. Discrete Contin. Dyn. Syst. Ser. B, to appear. 

[6] A. Baskaran, M. C. Marchetti, Nonequilibrium statistical mechanics of self-propelled 
hard rods, J. Stat. Mech. Theory Exp., (2010) P04019. 

[7] S. Bazazi, J. Buhl, J. J. Hale, M. L. Anstey, G. A. Sword, S. J. Simpson, I. D. Couzin, 
Collective Motion and Cannibalism in Locust Migratory Bands, Current Biology 18 

(2008) 735-739. 

[8] F. Berthelin, P. Degond, M. Delitala, M. Rascle, A model for the formation and 
evolution of traffic jams. Arch. Rat. Mech. Anal., 187 (2008) 185-220. 

[9] F. Berthelin, P. Degond, V. Le Blanc, S. Moutari, J. Royer, M. Rascle, A Traffic-Flow 
Model with Constraints for the Modeling of Traffic Jams, Math. Models Methods 
Appl. Sci., 18 Suppl. (2008) 1269-1298. 

[10] E. Bertin, M. Droz and G. Gregoire, Hydrodynamic equations for self-propelled par- 
ticles: microscopic derivation and stability analysis, J. Phys. A: Math. Theor., 42 

(2009) 445001. 

[11] F. BoUey, J. A. Cafiizo, J. A. Carrillo, Mean- field limit for the stochastic Vicsek 
model, Appl. Math. Lett., 25 (2011) 339-343. 

24 



[12] F. Brown, Micromagnetics, Wiley, New York, 1963. 

[13] E. Carlen, R. Chatelin, P. Degond, and B Wennberg, Kinetic hierarchy and propa- 
gation of chaos in biological swarm models, Phys. D, appeared online. 

[14] E. Carlen, P. Degond, and B Wennberg, Kinetic limits for pair-interaction driven 
master equations and biological swarm models. Math. Models Methods Appl. Sci., 
23 (2013) 1339-1376. 

[15] J. A. Carrillo, M. Fornasier, J. Rosado, G. Toscani, Asymptotic Flocking Dynamics 
for the kinetic Cucker-Smale model, SIAM J. Math. Anal, 42 (2010) 218-236. 

[16] H. Chate, F. Ginelli, G. Gregoire, F. Raynaud, Collective motion of self-propelled 
particles interacting without cohesion, Phys. Rev. E 77 (2008) 046113 (15 p.) 

[17] G. Q. Chen, C. D. Levermore, T. P. Liu, Hyperbolic conservation laws with stiff 
relaxation and entropy, Comm. Pure Appl. Math., 47 (1994) 787-830. 

[18] Y-L. Chuang, M. R. D'Orsogna, D. Marthaler, A. L. Bertozzi and L. S. Chayes, 
State transitions and the continuum limit for a 2D interacting, self-propelled particle 
system, Physica D, 232 (2007) 33-47. 

[19] I. D. Couzin, J. Krause, R. James, G. D. Ruxton and N. R. Franks, Collective 
Memory and Spatial Sorting in Animal Groups, J. theor. Biol., 218 (2002), 1-11. 

[20] F. Cucker, Er. Mordecki, Flocking in noisy environments, J. Math. Pures Appl., 89 
(2008) 278-296. 

[21] F. Cucker, S. Smale, Emergent behavior in flocks, IEEE Transactions on Automatic 
Control, 52 (2007) 852-862. 

[22] A. Czirok, E. Ben- Jacob, I. Cohen, T. Vicsek, Formation of complex bacterial colonies 
via self-generated vortices, Phys. Rev. E, 54 (1996) 1791-18091. 

[23] P. Degond, M. Delitala, Modelling and simulation of vehicular traffic jam formation, 
Kinet. Relat. Models, 1 (2008) 279-293. 

[24] P. Degond, A. Frouvelle, J-G. Liu, Macroscopic limits and phase transition in a 
system of self-propelled particles, J. Nonlinear Sci., appeared online. 

[25] P. Degond, A. Frouvelle, J.-G. Liu, A note on phase transitions for the Smoluchowski 
equation with dipolar potential, submitted. 

[26] P. Degond, A. Frouvelle, J-G. Liu, Phase transitions, hysteresis, and hyperbolicity 
for self-organized alignment dynamics, preprint. 

[27] P. Degond, J. Hua, Self- Organized Hydrodynamics with congestion and path forma- 
tion in crowds, J. Comput. Phys., 237 (2013) 299-319. 

[28] P. Degond, J. Hua, L. Navoret, Numerical simulations of the Euler system with 
congestion constraint, J. Comput. Phys., 230 (2011) 8057-8088. 

25 



[29] P. Degond, J-G. Liu, Hydrodynamics of self- alignment interactions with precession 
and derivation of the Landau-Lifschitz-Gilbert equation, Math. Models Methods 
Appl. Sci., 22 Suppl. 1 (2012) 1140001 (18 pages). 

[30] P. Degond, J-G. Liu, S. Motsch, V. Panferov, Hydrodynamic models of self-organized 
dynamics: derivation and existence theory. Methods Appl. Anal., to appear. 

[31] P. Degond, J.-G. Liu, C. Ringhofer, A Nash equilibrium macroscopic closure for 
kinetic models coupled with Mean-Field Games, submitted. arXiv:1212.6130. 

[32] P. Degond, S. Motsch, Continuum limit of self-driven particles with orientation in- 
teraction. Math. Models Methods Appl. Sci., 18 Suppl. (2008) 1193-1215. 

[33] P. Degond, S. Motsch, Large scale dynamics of the Persistent Turning Walker model 
offish behavior, J. Stat. Phys., 131 (2008) 989-1021. 

[34] P. Degond, S. Motsch, A macroscopic model for a system of swarming agents using 
curvature control, J. Stat. Phys., 143 (2011) 685-714 

[35] P. Degond, L. Navoret, R. Bon, D. Sanchez, Congestion in a macroscopic model of 
self-driven particles modeling gregariousness, J. Stat. Phys., 138 (2010) 85-125. 

[36] P. Degond, T. Yang, Diffusion in a continuum model of self-propelled particles with 
alignment interaction. Math. Models Methods Appl. Sci., 20 Suppl. (2010) 1459-1490. 

[37] M. L. Domeier, P. L. Colin, Tropical reef fish spawning aggregations: defined and 
reviewed. Bulletin of Marine Science, 60 (1997) 698-726. 

[38] A. Frouvelle, A continuum model for alignment of self-propelled particles with 
anisotropy and density-dependent parameters, Math. Mod. Meth. Appl. Sci., 22 
(2012) 1250011 (40 p.). 

[39] A. Frouvelle, J.-G. Liu, Dynamics in a kinetic model of oriented particles with phase 
transition, SIAM J. Math. Anal., 44 (2012) 791-826. 

[40] J. Gautrais, C. Jost, M. Soria, A. Campo, S. Motsch, R. Fournier, S. Blanco, G. 
Theraulaz, Analyzing fish movement as a persistent turning walker, J. Math. Biol., 
58 (2009) 429-445. 

[41] J. Gautrais, F. Ginelli, R. Fournier, S. Blanco, M. Soria, H. Chate, G. Theraulaz, 
Deciphering interactions in moving animal groups. Plos Comput. Biol., 8 (2012) 
el002678. 

[42] S. -Y. Ha, J.-G. Liu, A simple proof of the Cucker-Smale flocking dynamics and 
mean-field limit, Commun. Math. Sci., 7 (2009) 297-325. 

[43] S.-Y. Ha, E. Tadmor, From particle to kinetic and hydrodynamic descriptions of 
flocking. Kinetic and Related Models, 1 (2008) 415-435. 



26 



[44] E. P. Hsu, Stochastic Analysis on Manifolds, Graduate Series in Mathematics, Amer- 
ican Mathematical Society, 2002. 

[45] A. Khuong, G. Theraulaz, C. Jost, A. Perna, J. Gautrais, A computational model of 
ant nest morphogenesis, in "Advances in Artificial Life, ECAL 2011 - Synthesis and 
Simulation of Living Systems", MIT Press, 2011, pp. 404-411. 

[46] P.L. LeFloch. Entropy weak solutions to nonlinear hyperbolic systems under noncon- 
servative form, Comm. Partial Differential Equations, 13 (1988) 669-727. 

[47] A. Mogilner, L. Edelstein-Keshet, L. Bent and A. Spiros, Mutual interactions, po- 
tentials, and individual distance in a social aggregation, J. Math. Biol., 47 (2003) 
353-389. 

[48] J. Monod, Chance and Necessity: An Essay on the Natural Philosophy of Modern 
Biology, Alfred A. Knopf, New York, 1971. 

[49] S. Motsch, L. Navoret, Numerical simulations of a non- conservative hyperbolic system 
with geometric constraints describing swarming behavior, Multiscale Model. Simul., 
9 (2011) 1253-1275. 

[50] S. Motsch, E. Tadmor, A new model for self-organized dynamics and its flocking 
behavior, J. Stat. Phys., 144 (2011) 923-947. 

[51] V. V. Rusanov, Calculation of interaction of non-steady shock waves with obstacles, 
J. Comput. Math. Phys. USSR 1 (1961) 267-279. 

[52] J. Shen, Cucker-Smale flocking under hierarchical leadership, SIAM J. Appl. Math., 
58 (2007) 694-719. 

[53] J. Toner and Y. Tu, Flocks, Long-range order in a two-dimensional dynamical XY 
model: how birds fly together, Phys. Rev. Lett., 75 (1995) 4326-4329. 

[54] J. Toner, Y. Tu and S. Ramaswamy, Hydrodynamics and phases of flocks. Annals of 
Physics, 318 (2005) 170-244 

[55] Y. Tu, J. Toner and M. Ulm, Sound waves and the absence of Galilean invariance in 
flocks, Phys. Rev. Lett., 80 (1998) 4819-4822. 

[56] T. Vicsek, A. Czirok, E. Ben- Jacob, I. Cohen, O. Shochet, Novel type of phase 
transition in a system of self-driven particles, Phys. Rev. Lett., 75 (1995) 1226-1229. 

[57] T. Vicsek, A. Zafeiris, Collective motion, Phys. Rep., 517 (2012) 71-140. 



27 



